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Abstract 

We propose a simpfe equation for predicting self-diffusivity of fluids embedded in random ma- 
trices of identical, but dynamically frozen, particles (i.e., quenched-annealed systems). The only 
nontrivial input is the volume available to mobile particles, which also can be predicted for two 
common matrix types that reflect equilibrium and non-equilibrium fluid structures. The proposed 
equation can account for the large differences in mobility exhibited by quenched-annealed systems 
with indistinguishable static pair correlations, illustrating the key role that available volume plays 
in transport. 

PACS numbers: 66.10. Cb,46.65.+g,61.20.Ja 



Fluid transport in random porous media is central to a host of natural phenomena and 
technological applications, from the function of biological cells to the performance of mate- 
rials for membrane separations and heterogeneous catalysis. Although it is now appreciated 
that kinetic processes in these systems are intimately connected to their microstructures 
(pore volume, surface area, connectivity, etc.) determining the precise structure-property 
relationship for the transport property of interest remains a formidable challenge. 

One promising line of inquiry has been the exploration of simple model systems for which 
both transport properties and structure are amenable to theoretical analysis. In the case of 
diffusivity, earlier work has focused primarily on single particle transport through random 
configurations of static obstacles Q]. At low obstacle densities, particles show anomalous, 
sub-diffusive motion over intermediate time and length scales but recover normal diffusive 
behavior in the long time limit. On the other hand, sub-diffusive motion is observed for all 
times at high obstacle densities. Similar behavior is predicted for transport of single ions in 
disordered matrices of quenched charges 0|. 

Self-diffusivities of highly coupled fluids in random media pose additional challenges and 
are also of great interest. Very recently, insightful molecular dynamics simulation |4| and 
mode-coupling theory \^ studies have been presented for the hard-sphere (HS) fluid embed- 
ded in disordered matrices of obstacle spheres. This type of model, introduced originally by 
Madden and Glandt |(|, is referred to as a quenched-annealed (QA) system. To the best 
of our knowledge, no simple theoretical relationship between self-diffusivity of mobile parti- 
cles, densities of mobile and matrix particles, and matrix microstructure has been proposed 
for QA media. As a result, a general consensus on which physical parameters are most 
important for understanding mobility in these model porous materials is lacking. 

In this Rapid Communication, we introduce one such approximate relationship that we 
motivate by using a physical argument for how the matrix reduces the volume available for 
diffusion of the mobile particles. We test its predictions against the numerical results of 
molecular dynamics simulations for HS QA systems with matrices that reflect equilibrium 
and non-equilibrium structures, which serve as idealizations for two different classes of phys- 
ically realizable materials. We then demonstrate how the available volume in both types 
of matrices can be accurately predicted directly from equilibrium properties of the bulk HS 
fluid. Finally, we test our equation for diffusivity against simulation results of QA systems 
comprising Lennard- Jones (LJ) particles. For all models investigated, we find good agree- 
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ment with simulations over a wide range of parameters, with significant overpredictions of 
diffusivity occuring only for high matrix densities where limited connectivity of the available 
volume (not accounted for in our approximate equation) also hinders transport. 

The protocols that we use to generate static matrices represent two of the most commonly 
employed in studies of QA systems jj, 0, Q|. In the first, N M particles of diameter a are 
initially equilibrated in a volume V at temperature T. They are subsequently quenched (i.e., 
frozen in place) in an equilibrium configuration, and a "fluid" (F) of Np identical particles is 
added, equilibrated, and studied. This method of matrix generation is referred to as QA-M 
because the matrix (M) itself reflects an equilibrium configuration of density pm = Nm<7 3 /V 
for the bulk system. QA-M matrices serve as elementary models of amorphous solid materials 
prepared, e.g, by very rapidly cooling gel-forming suspensions of proteins or colloids 3]. 

The second matrix generation protocol we study is referred to as QA-FM because it 
involves first equilibrating Np + Nm identical particles of diameter a in a volume V at 
temperature T. Then, Nm of the particles are randomly selected and quenched to create the 
solid matrix. The other Np particles constitute the mobile "fluid" of density pp = Npa 3 /V. 
QA-FM matrices are idealized models for templated porous solids synthesized by depleting, a 
high density material of one of its components by, e.g., dissolution, reaction, or desorption 7|. 

Our approach for predicting the fluid self-diffusivities of these systems is motivated by 
three basic observations about their behaviors, (i) The total pair correlation functions of HS 
QA-M and QA-FM configurations, averaged over mobile and matrix particles together at 
the same total density pp + pm, are indistinguishable (ii) Yet, their fluid self-diffusivities 
generally differ significantly, even when compared at the same ratio of matrix to mobile 
particles. This has previously been interpreted as evidence that the differences in the dy- 
namics of QA media cannot be predicted based on static structural information alone j^. 
(iii) However, there is a key static property that distinguishes individual QA systems: the 
fraction of volume available to the mobile particle centers in the matrix, Vq/V pj. As we 
demonstrate below, differences in this quantity can largely account for the wide range of 
self-diffusivities exhibited by QA materials that are otherwise "structurally similar" . 

For example, contrast the behavior of a bulk HS fluid with density p = pp + pm to a 
HS QA system (produced by either -M or -FM protocols) with matrix and mobile particle 
densities of pm and p F , respectively. The binary collisions that mediate diffusion in these two 
systems are expected to be comparable, since their pair correlation functions and mobile- 
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particle thermal velocities are the same. Nonetheless, mobile particles in the QA system 
will diffuse over a shorter characteristic length scale per unit time compared to those of the 
bulk because the volume available for diffusion in the QA matrix, Vq, is smaller than V. If 
we assume that the characteristic length scales for diffusion are proportional to the cube 
root of the respective available volumes and, in turn, that the associated self-diffusivities 
are proportional to the square of the length scales, we arrive at the following approximate 
relationship: 

D(p F , p M ) « D(p F + p M , 0) x (V /V) 2/ \ (1) 

where D(p F , p M ) and D(p F + p M , 0) are the self-diffusivities of the mobile particles in QA and 
bulk fluid systems. As we show below, the appeal of this relation is its ability to predict the 
self-diffusivity of non-trivial QA systems from knowledge of well-characterized equilibrium 
properties of the bulk HS fluid. 

In HS QA-M systems, matrices are drawn from equilibrium configurations of the bulk HS 
fluid with density pm, and thus V /V is given by j| 

V /V = exp[-p c *(p M )/k B T}. (2) 

The excess chemical potential (relative to ideal gas) of the HS fluid, /i ex , can be obtained an- 
alytically from, e.g., either scaled particle theory |S| or the Carnahan-Starling (CS) equation 
of state |9J] (here, we adopt the latter, which is very accurate for pu < 0.9). 

In HS QA-FM media, however, the matrix structure is different than that of an equilib- 
rium HS fluid because it is quenched in the presence of mobile particles. Nonetheless, as Van 



Tassel and co-workers 



have demonstrated, because of the specific protocol by which 
"non-equilibrium" QA-FM matrices are created, one can still apply equilibrium liquid-state 
approaches (e.g., integral equation theories) to estimate Vq/V. Here, we introduce an accu- 
rate information-theory (IT) based stategy for accomplishing this task. As we show, its main 
advantage is that it only requires knowledge of pm, Pf, and the pair correlation function g(r) 
of the bulk HS fluid. 

Our strategy implicitly uses the fact that the pair correlations of the HS QA-FM system 
are indistinguishable from those of an equilibrium HS fluid with density p F + pM- Since the 
QA-FM particle identities (matrix or mobile) are randomly distributed, one can write down 
an exact expression for Vq/V in the QA-FM matrix in terms of equilibrium properties of 



4 



bulk HS fluid: 

Pm 



Vo/V = J2^i(pM + Pf) x 

i=0 



1 + 
Pf 



(3) 



Here Ilj(pM + pp) represents the probability that a randomly placed spherical window (with 
radius equal to one particle diameter) in a bulk HS fluid of density pf + Pm will contain 
precisely i particle centers (see Fig. |TJ). The quantitiy [1 + Pm/pf] 1 is the probability that, 
in the structurally equivalent QA-FM system, all of the i centers in the window would be 
fluid particles. 

IT provides expressions for the n^, 

= exp[A + Aii + \ 2 i 2 ]/i\, (4) 

which maximize a relative information entropy subject to some experimental constraints 
In particluar, Ao, Ai, and A2 are Lagrange multipliers determined (see, e.g., [ill) by imposing 
the normalization condition J2i n« = 1 and the first two moments of the window occupancy, 
i and i 2 , 

i = 4np/3, 



1 



2 



1 + p 2 dr dv'g{\r-v'\). (5) 



Here, p = Na 3 /V, and the subscript w in the last expression indicates that the integrals 
are constrained to the spherical observation window. Fig. |21 shows a numerical comparison 
of Vo/V predicted by Eq. |3] and |U to the "exact" results for Vq/V obtained from applying 



the available space algorithm of Sastry et al. 12j to HS QA-FM matrices. As is evident, the 



simple IT approach provides accurate predictions over a wide range of matrix parameters. 

Having established means for estimating Vo/V of both QA-M and QA-FM matrices, we 
now only require an expression for the self-diffusivity of the bulk HS fluid, D(p, 0), in order to 
predict HS QA diffusivity using Eq. [T] For this, we adopt Speedy's empirical fit to molecular 
dynamics simulation data 

D(p, 0) = 4 ( 1 - -n^] [l + P 2 (0.4 - 0.83p 2 )l , (6) 



p \ 1.09 

where A = 3/(8-\/7r). We implicitly non-dimensionalized D(p, 0) of Eq. and all other self- 



diffusivities in this study, by a /r, where r = aym/k^T, k-Q is the Boltzmann's constant, 
and m is particle mass. 
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To test the validity of Eq. ^ we performed molecular dynamics simulations for HS 
QA-M and QA-FM systems in the microcanonical ensemble using N F = 1500 mobile particles 
and a periodically-replicated simulation cell. Runs were performed at different values of Nm 
and V/a 3 to obtain results for specific combinations of p-p and pM- We extracted self- 
diffusivity of the mobile particles D(p F , p M ) by fitting the long-time (t ^> 1) behavior of the 
average mean-squared displacement to the Einstein relation (Ar 2 ) = 6Dt. 

Fig. 01 provides a fairly comprehensive comparison of D(pp,pM) from Eq. ^to the results 
of our molecular dynamics simulations. Interestingly, the predictions show semi-quantitative 
agreement with simulations of QA systems for matrix densities in the wide range < pu < 
0.25, which includes systems where matrix particles exclude mobile particle centers from over 
70% of the total volume (see, e.g., Fig. |2J). Given that matrix particles reduce D(p F ,p M ) by 
more than an order of magnitude across this range of densities, the success of Eq. ^ argues 
that available volume plays a primary role in controlling the single-particle dynamics. For 
matrix densities greater than pM = 0.25 (i.e., systems with less than 30% of volume available 
to mobile particle centers), Eq. [T] captures the qualitative trends, but it systematically 
overpredicts the simulated self-diffusivities. This overprediction is, of course, expected, 
given that Eq. ^ only accounts for the reduction of available volume and not the fact that 
available volume also becomes highly disconnected (i.e., some pockets of available volume 
are inaccessible or are accessible by only small number of paths) at high matrix densities Q , 
which acts to further hinder transport. However, for systems of physical relevance (porous 
catalysts, membranes for separations, transport gels for drug delivery, etc.) which rely on 
rapid diffusion for their functionality, the precise connectivity of the available volume will 
be less important. 

In order to test whether the connection between diffusivity and available volume holds 
more generally for fluids confined to random media, we also performed molecular dynam- 
ics simulations of QA-FM systems comprising LJ particles (truncated and shifted with a 
quadratic function in r to insure that both the potential and its gradient vanish at 
r cut = 2.5). The details of the simulations are identical to those of the HS QA systems 
described earlier, except that the equations of motion were integrated via the velocity Ver- 
let algorithm [1^, and N-p = 1000 mobile particles were considered. The goal was to test 
if one could employ techniques commonly used in thermodynamic perturbation theory to 
map the QA L J system onto an equivalent QA HS system, and then use Eq. [T] to predict 
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the self-diffusivity. The specific mapping that we used in this study is a Boltzmann fac- 



tor criterion 
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17| . which determines the temperature-dependent "effective" HS diameter 
cr(T) (and hence the corresponding reduced matrix and fluid denisties, pm and p F ) of the LJ 
system through the following relationship u (r = a) = k^T, where Uq is the repulsive part 
of the Weeks- Chandler- Andersen decomposition of the pair potential The hypothesis 
is that LJ QA systems will exhibit similar reduced diffusivities D(p F , pu) as HS QA systems 
when compared at the same values of p F and pM- 

As a first test of this idea, we compare in Fig. the diffusivities D(p F ,p M ) of the HS 
QA-FM and LJ QA-FM systems obtained by our molecular dynamics simulations. In this 
plot, the LJ QA-FM system is at a temperature T = e/k-Q, where e is the characteristic 
energy of the LJ pair potential Except for at very low fluid densities (conditions for 
which the LJ fluid structure is not accurately approximated by an equivalent HS reference 

rn 

fluid (l9|)> the mapping brings the dynamics of the LJ and HS QA systems into excellent 
agreement. 

A final test is to check whether the simple aforementioned mapping allows direct predic- 
tion of L J QA-FM dynamics at other temperatures using Eq. ^ Fig. |3Jd shows the compar- 
ison of the predicted D(p F ,p M ) versus the results of our molecular dynamics simulations at 
T = 3e/k-B- The agreement between the predictions and the simulations again confirm the 
pivotal role that available volume plays in controlling the single-particle dynamics of fluids 
in porous media. 
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FIG. 1: Schematic of a QA system. Equi-sized filled and empty circles are matrix and fluid 
particles, respectively. The bigger circle is a window with radius equal to the one particle diameter. 
Hi is the probability that it will be occupied by exactly i matrix or fluid particle centers. 
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FIG. 2: The fractional available volume Vq/V in HS QA-FM matrices computed using the IT 
approach of Eq. and (curves) and "exact" results (filled circles) using the available space 
algorithm of Ref |l2| • 
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FIG. 3: Fluid self-diffusivity D versus fluid density p F for (a) HS QA-FM and (b) HS QA-M 
systems with matrix densities pm = 0.0, 0.05, 0.15, 0.25, 0.30, and 0.35. Curves are predictions of 
Eq. ^and circles are molecular dynamics simulation results. 
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FIG. 4: (a) Fluid self-diffusivity D versus fluid density p F for HS QA-FM (circles) and LJ QA- 
FM (diamonds) systems obtained via molecular dynamics simulations. The LJ QA-FM system 
is at T = e/k-Q. Matrix densities of pm = 0.05, 0.15, 0.25, and 0.35 are presented, (b) Results 
for the LJ QA-FM system at T = 3e/&B : molecular dynamics simulations (circles) and Eq. ^ 
(curves). Matrix densities of pm = 0.043, 0.128, 0.214, and 0.3 are presented. In both (a) and 
(b), matrix and mobile densities for the LJ QA-FM systems are defined as pm = AmO"(T) 3 /V and 
Pf = Npa(T) s /V , respectively. The effective HS diameter cr(T) is determined by a Boltzmann 
factor criterion described in the text. 
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